#! /bin/sh
WIDTH=600
HEIGHT=400
WIDTHOFF1=0
WIDTHOFF2=640
HEIGHTOFF1=50
HEIGHTOFF2=530

deconparams=deconparams3
fold=8
space=3		# 3 null traces between panels
tpow=1		# gain power before examining autocorrelations
		# tpow=2 would often be better for field data
mix=1,1,1
pnoise=.0005

###PURPOSE###
# We use the same data described in the first demo.  But here
# we illustrate how to set up a series of tests to determine the
# optimal gap for prediction error filtering.


# First make the synthetic data for the deconvolution demo.
######START OF MODEL######
I=${CWPROOT}/include
L=${CWPROOT}/lib

make

./traces |
suaddhead ns=512 |
sushw key=dt a=4000 |
suaddnoise sn=2000 |
suconv filter=1152,-384,-904,288,174,-34,-12 >modeldata

rm traces
######END OF MODEL######

# Plot the model 
suxwigb <modeldata title="Data: 64ms reverbs" \
	windowtitle="Data" \
	label1="Time (sec)" label2="Trace" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF1 ybox=$HEIGHTOFF1 &
 
# Plot the autocorrelation with gain
sugain <modeldata tpow=$tpow |
suacor nt=31 |
suxwigb title="Autocorrelation" \
	windowtitle="AutoCorr" \
	label1="Time (sec)" label2="Trace" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF2 ybox=$HEIGHTOFF1 &


#!# Determine range of decon parameters (use milliseconds to run loop)
gapmin=40 gapmax=40
lenmin=100 lenmax=180
gapincrement=0 gapexpand=1
lenincrement=10 lenexpand=1
tmin=0.5 tmax=2.0	# auto-corr window
f1=8 f2=12 f3=95 f4=100	# bandpass frequency corners

panelnumber=0  # test panel number


### Determine ns, and dt from data (for sunull)
nt=`sugethw ns <modeldata | sed 1q | sed 's/.*ns=//'`
dt=`sugethw dt <modeldata | sed 1q | sed 's/.*dt=//'`

### Convert dt to seconds from header value in microseconds
dt=`bc -l <<-END
	scale=3
	$dt / 1000000
END`

# Loop over deconvolution parameters
>deconpanel  # zero output files

# zeroth panel is w/o decon
cp modeldata deconpanel
sunull ntr=$space nt=$nt dt=$dt >>deconpanel

echo "Making deconvolution test panels.  Takes a few seconds ..."
echo "Zeroth panel is data without deconvolution."

echo "Test minlag maxlag tmin tmax" >$deconparams
gap=$gapmin
len=$lenmin
while [ $gap -le $gapmax -a $len -le $lenmax ]
do
	panelnumber=`expr $panelnumber + 1`
	minlag=`bc -l <<-END
		scale=5; $gap / 1000
	END`
	maxlag=`bc -l <<-END
		scale=5; $minlag + $len / 1000
	END`

	sugain <modeldata tpow=$tpow |
	supef pnoise=$pnoise mix=$mix minlag=$minlag maxlag=$maxlag mincorr=$tmin maxcorr=$tmax |
	sufilter f=$f1,$f2,$f3,$f4  >>deconpanel
	sunull ntr=$space nt=$nt dt=$dt >>deconpanel

	echo "$panelnumber $minlag $maxlag $tmin $tmax" >>$deconparams
	gap=`bc <<-END
		scale=5
		$gap*$gapexpand + $gapincrement
	END`
	len=`bc -l <<-END
		scale=5
		$len*$lenexpand + $lenincrement
	END`
done

echo "The parameter values are recorded in $deconparams:"
cat $deconparams

### Plot deconvolution test panels
f2=0
d2=`bc -l <<-END
	scale=5
	1/($fold + $space)
END`

### Additional display gaining and plot
### Some might like tpow=1 pbal=1 instead of gpow=.5, etc.
sugain <deconpanel gpow=.5 |
suxwigb f2=$f2 d2=$d2 title="PEF Deconvolution Test" \
	windowtitle="DeconTest" \
	label1="Time (s)"  label2="Deconvolution Test Number" \
	f2num=0 d2num=1 \
	wbox=1200 hbox=$HEIGHT xbox=25 ybox=$HEIGHTOFF2 &

exit 0
